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We give a hydrodynamical explanation for the chaotic behaviour of a dripping faucet 
using the results of the stability analysis of a static pendant drop and a proper orthogonal 
decomposition (POD) of the complete dynamics. We find that the only relevant modes are 
the two classical normal forms associated with a Saddle- Node- Andronov bifurcation and a 
Shilnikov homoclinic bifurcation. This allows us to construct a hierarchy of reduced order 
models including maps and ordinary differential equations which are able to qualitatively 
explain prior experiments and numerical simulations of the governing partial differential 
equations and provide an explanation for the complexity in dripping. We also provide a 
new mechanical analogue for the dripping faucet and a simple rationale for the transition 
from dripping to jetting modes in the flow from a faucet. 



1. Introduction 

Almost since the beginning of the modern revolution in nonlinear dynamics and chaos, 
the dripping faucet has served as a paradigm of chaotic dynamics l|Shaw 1984|l . To explain 
the transition to chaos, various ad-hoc mechanical models based on variable-mass and 
spring systems have been used to derive return maps and Poincare sections for the time 
between droplet emissions. Nearly twenty years on, this empirical approach still continues 
(see ( |Kiyono et al 1999) for a recent example along with a review of earlier work) , while 
the connection to hydrodynamics remains tenuous. In contrast, over the last decade, the 
hydrodynamical approach has been used with great success in studying the pinch-off 
of a single drop using a combination of theory, numerical simulation and experiment 
(for a review see ( |Eggers 1997| l). Thus it is natural to bring this understanding of the 
hydrodynamics of drop emission to bear upon the complex dynamics of dripping and 
understand its origins. Two recent attempts in this direction are Fuchikami ct al. 199^ 
and Am bravaneswaran et al 20001 who simulate the chaotic dynamics of dripping using 
simplified lubrication-like approximations of the Navier-Stokes equations. However, the 
mechanistic origins of the complexity in dripping still remain elusive. 

Our goal in this paper is not as much to match the images of chaotic dripping seen 
in experiments or numerical simulations as to see how we might qualitatively explain 
what seems to be a complex dynamical process as simply as possible. With this in mind, 
we derive a series of models culminating in a mechanical analogue, all of which cap- 
ture the salient features of periodic and chaotic dripping and jetting, thus bringing a 
focus on the physical mechanisms responsible while providing a guide to further work on 
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related problems. Some preliminary results of this nature were first announced in a con- 
ference proceedings l|Coullet et al 20d(Hl a few years ago. In §2, we revisit the equations 
of equilibrium for a static pendant drop. In §3, we describe a Lagrangian model based on 
lubrication theory for the evolution of the drop and use it to carry out a linearized anal- 
ysis of the static drop shapes as a function of their volume. In addition, we also present 
solutions of the complete dynamical evolution associated with dripping for a range of flow 
rates and evaluate the bifurcation diagram that shows the transition to chaos. In §4, we 
analyze these chaotic solutions using proper orthogonal decomposition and a time-delay 
method, both of which suggest that a low dimensional model should suffice to explain the 
complexity of dripping. In particular, we show that the dynamics can be recast in terms 
of a hierarchy of low-dimensional models, i.e. coupled ordinary differential equations and 
maps, in §5 which are qualitatively consistent with prior experiments and yield a simple 
physically consistent picture of the complex dynamics in terms of well known concepts 
in dynamical systems. Guided by these models, we propose a new mechanical analogue 
of the dripping faucet in §6. We conclude the paper with a brief discussion in §7 of how 
our models can rationalize the transition from dripping to jetting. 



2. Equations of equilibrium 

We start with the case where the flow rate is very small, so that a drop remains 
attached to the faucet until its volume exceeds a threshold Vc- For a faucet of radius 
R that is sufficiently small, drops with a volume V < Vc are stable and axisymmet- 
ric ( |Padday and Pitt~1973j which case we will restrict ourselves to. The shape of such 
a pendant drop is determined by minimising the sum of its gravitational and surface 
energy subject to the constraint of constant volume. Alternatively, we can write the 
dimensionless equation for the balance of forces normal to the interface as 

dO cos 9 

= -z, (2.1) 

as r 

where the drop interface {r{s), z{s)) is determined by the kinematic conditions 

dz 



— cos 9 

, (2-2) 
— = sm 
ds 

Here the variables r, s, 9 and z are defined in Fig. QJi). The characteristic scales 
are : = \/^ IqP for length, mo = pl% for mass, = {V j pg^Y^'^ foi" time, where F, 
p, g are the surface tension, density of the fluid and gravity respectively. For water at 
2Q°C, Zo = 0.27cTO,TOo = 0.020g.to = 0.017s. The boundary conditions at the bottom 
of the drop are r(0) = 0, 0(0) = 7r/2 and z(0) = Pb/ pg, where Pb is the unknown 
hydrostatic pressure at this position and serves as a control parameter that describes 
the family of stationary drops. Choosing a value for Pb, we integrate (|2. 112.2(1 as an 
initial value problem until r = R [R = est) and find the corresponding drop volume 
V = J nr^dz. For a given Pt, when R < 1 the condition r = R may be satisfied for up to 
three different drop volumes and lengths URiera fc Risler 2002|l . These pendant drops are 
shown in Fig. for different faucet radii R. However, the only stable stationary drop 
corresponds to the branch starting at the origin and ending at the first turning point 
( |Padday and Pitt 1973| ) which denotes the critical drop volumes Vc for which the weight 
of the drop is just balanced by the capillary forces. The corresponding shapes have a 
waist where the dynamic instability associatedwith eventual pinch off first starts. 
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Figure 1. a) Schematic of the static drop, b) The dimensionless hydrostatic pressure at the 
bottom of the drop Pb as a function of the dimensionless volume V of the drop for different values 
of the dimensionless faucet radius R = 0.5 (left curve) and R = 1.0 (right curve). Solutions of 
12. H and 12.21 1 yield the P — V curves shown, with the the solid line corresponding to stable 
static drops and the dashed line corresponding to unstable stationary drops, 



3. Equations of motion 

To understand the mechanism of this hnear instabihty dynamically, we consider the hy- 
drodynamical equations linearized about a stationary solution. Instead of using the com- 
plete Navier-Stokes equations for Eulerian one-dimensional lubrication theories HEggers fc Dupont 1994| 
we simplify the analysis by using a new lubrication model embodied in a Lagrangian ap- 
proach for the fluid (^Fuchikami et al. 1999|l . The inherent assumptions in this are the 
following: (a) the drop remains axisymmetric during its motion, (b) the radial compo- 
nent of the fluid velocity is negligible compared to the axial component which depends 
only on z, (c) there is no overturning of the interface r(z). These assumptions are asymp- 
totically valid for slender drops of large viscosity, but recent simulations of the resulting 
low-order equations l|Fuchikami et al. 1999|l have shown good agreements with experi- 
ments even for fairly squat drops of low viscosity. The above assumptions lead to the 
conclusion that there is no exchange of fluid between neighboring horizontal slices of the 
drop, so that the volume of a slice is constant during the motion and can be treated as 
a Lagrangian variable. This yields a formulation that is essentially equivalent to earlier 
Eulerian description. Explicitly, the volume between the planes Zb{t) and Za{t) is 

aza,Zb,,t)= f '^\r{^,tfdC (3.1) 

where r is the radius of the drop. In terms of the Lagrangian variable ^(za, z^, t), we 
can write the kinetic energy Ekim potential energy Ug and surface tension energy Ur of 
the system as 
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Figure 2. The evolution of the drop as a function of time computed by solving 13.61 . The 
snapshots are separated by 0.5 unit of normalised time, unless marked otherwise. Parameters 
values for the simulation are 7? = 1, ?7 = 0.02, g = 1, 7 = 1, hq = 0.1. 



^kvn — 2 Jo \ dt ) "''^ 

U, = -pgjt^'^zitm (3-2) 
Ur = r ^w + L^d^ 

Here £,o{t) is the total volume of the drop at time t, and a prime corresponds to a partial 
derivative with respect to the Lagrangian variable ^. Then we can write the Lagrangian 
of the system as 

C = Ekra -Ug-Ur (3.3) 
The effect of viscosity is then expressed using the Rayleigh dissipation function 







Then, Lagrange's equation for the system is 

d dC dC 1 dEk,, 



dtdv dz ' 2 dv ^ ^ ^'^^ 

For the purposes of computation we follow the method of'Fuchika mi et ai. 19991 and 
discretize the Lagrangian spatially so that the drop is effectively sliced into N disks; each 
disks is characterised by 3 variables : the position z.^, the velocity Vi = ^ and the mass 
rrii . A disk is divided in two when its relative thickness (ratio of its height to its width) 
exceeds a threshold (0.05 in a typical simulation) to keep it slender. This occurs, either 
because of the added volume at the faucet or the change in the shape of the drop in 
the necking area. Two consecutive disks are combined when their relative thickness falls 
below a different threshold (0.075). The splitting and merging are done so as to conserve 
volume and momcntumf . This auto- adaptive mesh leads to a variation of N between 50 
and 4000 during a typical simulation (see Fig. El and Fig. I^J. 

The discretized Lagrangian then yields the TV equations of motion 

f Other quantities are not conserved as accurately, as in IFuchikami et al. 19991 
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Figure 3. Time interval between drops as a function of the exit velocity vo, obtained by solving 
H8.6^ . with R = 0.5,7? = 0.02, jr = 1,7 = 1. The qualitative features that are apparent are the 
oscillation of the time between drops as a function of the exit velocity vo, the single period-2 
oscillations before the chaotic region and finally the transition to chaos via period doubling and 
via a boundary crisis. All these features are robust functions of parameters, i.e. they persist for 
a range of values. In the inset (i) , we show the bifurcation diagram associated with increasing 
flow rate and in (ii) the corresponding diagram with decreasing flow rate; we observe that there 
is hysteresis. 



d dC dC 1 dE 

kin - -1 AT / n n\ 

dt dvi dzi 2 dvi ' ' 

These equations of motion are then integrated using a fifth-order Runge-Kutta method 
with an adaptive timestep that is based on the local truncation error (see Fig. I^J . The 
initial condition is either a stable stationary solution of the equation of equilibrium 
H2.1I2.2|I or a hemisphere. 

To determine the stability of the static solution we use the solutions determined in 
§2 via a shooting method (since there are no analytic solutions for the static shape) 
and substitute them into the linearized equation of motion H3.6|) in the neighbourhood 
of the stationary solutions of 12.1I2.2|I and determine the eigenvalues uji,i = 1,N of 
the resulting system. The value of the real part and imaginary part of the 3 largest 
eigenvalues are plotted in Fig. When V < Vc, 5R[wi] < 0, so that these drops are 
stable. As ^ Vc two complex conjugate eigenvalues become real by colliding on the 
real axis (Fig. (|SJ)) corresponding to the critically damped oscillator. As V increases 
further, the eigenvalues split; but remain on the real axis. One of the eigenvalues then 
moves away from the imaginary axis and the other moves towards it eventually reaching 
the origin when V — Vc- This is consistent with Fig. where we show the volume of 
the drop as a function of P},. The collision of two branches of stationary solutions one 
stable and one unstable and the disappearance of the stationary solution when V > Vc 
is a characteristic of a saddle-node bifurcation. 

At low flow rates, the shape of the drop is always close to that of a stationary drop. 
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Figure 4. Evolution of the real and imaginary part of the largest three eigenvalues 
(SR[a;i], ^[tJi]) as a function of the scaled drop volume V shows that 5R[lji] ~ 0.15R[tJ2] 
(7? = 0.5, = 0.02, = 1,7 = 1). 
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Figure 5. Sketch of the location of the first few eigenvalues when V ^ Vc- The first "nose" of 
the Pi, — V curve corresponds to the collision between the unstable "saddle" and stable "node" 
solutions and leads to a saddle-node bifurcation. 



The spectrum of the linearized stable solutions with V ^ Vc, as shown in Fig. I^J in- 
dicates that only two modes are dynamically relevant because the others are rapidly 
damped. These modes correspond to a) an oscillatory damped mode given by the largest 
two complex conjugate eigenpairs and b) a Saddle-Node bifurcation (Fig. ((SJ). To under- 
stand the dynamics of drop emission and relaxation, we now turn towards the numerical 
solution of the equations of motion (|3.6I) with the eventual goal of justifying a low order 
approximation for this complex process. 



4. Galerkin proper orthogonal decomposition 

Proper Orthogonal Decomposition (POD) ( [Holmes et al 1998) \ allows us to project 
a a high/infinite dimensional system onto a finite number of basis functions or modes. 
This may be done using the data from the direct numerical simulations, time series etc. 
and allows us to try and understand the qualitatively important features of a complex 
process. The optimality of POD derives from the fact that a truncated POD describes 
typical members of the ensemble better than any other decomposition. We use POD 
to determine if a low mode approximation of the dynamics of dripping is valid. As an 
example we carry out a POD in a regime with period-2 dripping {vq = 1.025 <C Io/tq-, flow 




Figure 6. A proper-orthogonal decomposition (POD) of the numerical simulation of the hy- 
drodynamical equations 13.611 for the dripping faucet {R = 0.5, rj = 0.02 and vq = 1.025). a) 
Singular values are plotted for each "mode", b) The first three modes from the POD superposed 
on the average drop shape: 1 corresponds to Vi+ra, 2 corresponds to V2+rav and 3 corresponds 
to V3 + ra(the dashed line corresponds to the interpolated average drop shape Va) computed as 
described in the text. Here Vi corresponds to the ith column of the matrix POD as indicated in 
the text. We see that just the first three modes capture much of the dynamical behavior of the 
faucet. 

rate q = ttR^vq and R — 0.5) where there is noticeable difference in the time required for 
the formation of large and small drops, allowing us to clearly distinguish the two periods. 

The shape of the drop r(z, t) is recorded periodically and a regular spaced mesh based 
on the maximum drop length is used to interpolate the typically irregular mesh. If a 
point in the interpolated regular mesh is outside the drop we set the radius r = 0. This 
mesh is used to compute both the average shape of the drop and the variation from 
this shape. This procedure allows us to follow the variance of the shape and leads to 
a rectangular matrix M of size m x n corresponding to m snapshots of the variation 
of the drop shape on a regular spaced mesh containing n points. The POD of M is a 
decomposition of the form M — UDV'^, where is an n x n matrix and U an m x n 
matrix. The columns of the matrix V form an orthonormal basis corresponding to the 
spatial modes of the drop. The columns of the matrix U also from an orthonormal basis 
and correspond to the normalised projection of each snapshot onto the spatial modes 
of the POD, i.e. MV — UD where the diagonal matrix D corresponds to the singular 
values. 

In Fig. |5| we show that the singular values of the first two modes are distinct from 
the other clustered modes for the particular choice of parameters indicated earlier. This 
separation of the singular values is strongly dependent on the scaled viscosity 77; as 77 
decreases all the singular values cluster together. The first mode ( 1 in Fig.E)) is associated 
with the process of drop growth and emission (Fig. [T] solid line). The second mode (2 in 
Fig. O corresponds to the largest damped oscillatory mode (Fig. [7| short dashed line). 
Since the higher modes are rapidly damped, we see that the dynamics of the dripping 
faucet can be well approximated using only the first two modes of the basis obtained 
using POD in this period-2 regime. Similar results are obtained in different regimes. 

We note that a priori there is no direct connection between the linearized modes and 
the POD modes. Here we distinguish the two uses of POD. Traditionally, it is used 
a method of simplifying the complex dynamics of high dimensional system. Here our 
goal is slightly different; what we want to show is that a low dimensional description is 
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Figure 7. Projection of the dynamics embodied in the solution to 13.61 . onto the first three 
POD modes. The solid line corresponds to the first mode UiDi, short dashed line the second 
U2D2 and long dashed and dotted line the third U3D3, with Uj being the jth column and Dj 
being the diagonal entry in the matrix POD as indicated in the text. 

possible, and in following we use the linearly stable modes to construct a hierarchy of 
low dimensional models. 

We can also reconstructed the phase space of the dynamics by using the time delay 
method (^Eckmann and Ruellc 1985). We use the radius of the drop at a given position 
above the usual location of the pinch-off so that the variable remains continuous during 
the process. We choose different delays Tq = 0, Ti, ...,Tjv-i to define our new variables 
fk — f{i + Tfc) which gives us a A'^-dimensional system. The minimum dimension N of 
the system is chosen so that the flow doesn't cross itself in the phase-space and the 
delays are chosen somewhat arbitrarily in order to obtain a reasonable projection. We 
find that = 3 is sufficient to capture the dynamics suggested by both the linear analysis 
and the POD. In the phase-space generated by the time-delay method shown in (Fig. 
(jHJ we see two qualitatively different regions: a large excursion corresponding to the 
dynamics that lead to drop pinch-off, and a much more compact region corresponding to 
the damped oscillations following the pinch-off event, that eventually leads the orbit to 
the neighbourhood of the saddle-node area whence it escapes again. 

5. Low dimensional models 

We now consider the dynamics of dripping using simplified models. The time scale 
for the formation of a pendant drop is ~ Once the volume of the pendant drop 
is Vc, it becomes unstable and pinches off in finite time. This process occurs in a time 
Tn ~ R^JpRlT (Clanct and Lasheras 1999) which is much shorter than the time for a 
new drop to form, and independent of the flow rate. Therefore the dynamics of dripping 
will not be affected by the details of drop pinch-off. 

After pinch-off, the remaining liquid recoils due to capillary forces, oscillating with a 
characteristic frequency / ^ ^JtJJjjV) which varies with the volume V and shape of the 
residual drop. For a given flow rate the volume of the pendant drop grows steadily, and 
this frequency gradually decreases even as the oscillations are damped out by viscous fluid 
motions at a rate l/rd ~ yJvJ^^^ (Fig. EJ. For very small flow rates, these oscillations 
are completely damped out by the time the pendant drop attains the critical volume 
Vc, so that in this case, droplets are emitted with a constant periodicity. As the flow 
rate is increased, these partially damped oscillations modify the onset of the instability 
via a saddle-node bifurcation. Equivalently, the dimensionless ratio of the flUing time 
to the damping time Tf/Td advances or delays the onset of necking and is responsible 
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Figure 8. Reconstruction of the flow by the time delay method obtained by solving 13. 6t 
numerically with parameter values corresponding to those for Fig |21 The radius of the drop 
r{t + it) = r{0.5,t + iT),i = 1,2,3 so that it remains continuous during the process, and the 
time delay r = 1.0. We observe a long excursion followed by a damped oscillations before the 
orbit returns to the neighbourhood of the saddle- node point. The arrow indicates the direction 
of the flow. 

for the variation of the periodicity (or lack thereof) of drop emission. For example, as 
the flow rate is gradually increased, the constant periodicity "drop-drop" gives way to 
a "drop-drip" scenario via a period-doubling bifurcation as follows. Once the pendant 
drop reaches the critical volume T4, a large droplet "drops" leading to a highly elongated 
residual filament. If the flow rate is large enough so that the oscillations are not completely 
damped out, the next droplet will start the necking process when V < Vc oi V > Vc 
depending of the direction of the oscillation close to the critical volume, so that a smaller 
or larger droplet "drips" . This leads to a smaller or larger residual drop whose oscillations 
will be damped out much sooner or later, thereby (possibly) allowing the pendant drop 
reach its maximum size before or after it "drops" , and so on. The temporal spacing 
between two drops arises from the damped oscillations close to the critical volume which 
enables the drop to become smaller or larger before the pinching off regime. 

Based on these observations, and the linear stability analysis consistent with the POD 
we build a simple low dimensional model. Let X be a state variable characterising the 
drop in a potential determined by the gravitational and surface energies as shown in 
Fig. 1^1 Giving a physical interpretation of X is tricky; it can be related to the center 
of gravity of the drop, as used for example by [Kiyono et al 2003| but only indirectly 
because different shapes can have the same center of gravity. A better interpretation 
of X relates it to the projection of the complete shape on a spatial mode of the drop 
and indicates the location of the drop relative to its stationary solutions as shown in 
Fig. 1^1 in the spirit of qualitative behaviour of dynamical systems. In the vicinity of the 
minimum of the potential energy the drop oscillates stably with a frequency determined 
by the curvature of the potential and damping rate determined by the spectrum of the 
linearized equation of motion. Close to the saddle point S the drop can oscillate stably 
or start necking, depending on the size of the perturbation. Indeed the height of the 
barrier represented by the unstable branch U is the energy for drop nucleation. A simple 
model potential characterising this \s E = a{y)X + ^^^^X^ , where V is the drop volume. 
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Figure 9. Potential corresponding to the capillary-gravitational energy E — a{V)X + b{V)X^ /3 
of the stable stationary pendant drop as a function of the drop volume V and the state of the 
drop X. The stable position (solid line) for the pendant drop corresponds to the local minimum 
of the potential, while the dashed line denotes the unstable solution. 

This potential evolves with the drop volume so that the stable solution disappears via 
a Saddle-Node bifurcation when a{Vc) — 0. Using these facts, we may then describe the 
dynamical evolution of the drop via the equations 



The first equation describes the dynamics of the state variable X{t) with a scaled 
damping r]{V) in the potential E, while the second equation quantifies the increase in 
the drop volume due to a flow rate q. The first equation of (|5.1I) charaterizes a simple 
damped oscillator; the nonlinearity is the normal form for the saddle-node bifurcation 
while the damped oscillations arise when we account for the stable oscillations of the 
recoiling drop. Thus the equation is a natural consequence of the stability analysis and 
the POD of the governing PDE. A mechanical analogue of H5.1|l is a damped particle 
moving on a curved surface forced in the direction V with a velocity q but free in the 
orthogonal direction X. To ensure that the potential characterizes the drop close to its 
stationary solution the spectrum of our dynamical system H5.1(l must approximate that 
obtained from the linear stability analysis of the stationary pendant drop described in 
§3. This leads to the following approximations for the functional form a, b, rj for V ^ Vc : 



where a, (3, 7 and 5 are fitting parameters computed to match the spectrum of the 
H3.6(l . The fitting functional form for (|5.2|l of a, 5, 77 are chosen to be simple and to have 
a minimum of parameters. However other functional forms do not change the qualitative 
behaviour of the solution to H5.1|l . The state variable X eventually diverges when it goes 
over the unstable branch U or when V > Vc which event corresponds to the ejection 
of a drop. The system is then reset bringing back X and V close to a stable solution 
in the same way as the drop goes back to a volume V < Vc and close to a stable drop 



{ 



dttX + 7j{v)dtX 

dtV 



a{V)+b{V)X^ 

q 



(5.1) 




(5.2) 
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Figure 10. (a) Total volume at the break up time as a function of the exit velocity no and 
(b) the ratio of the volume of the emitted drop to the total volume of the drop. To a first 
approximation, this ratio varies linearly with the exit velocity vo- 
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Figure 11. Bifurcation diagram derived from (a) the full hydrodynamic model 13. 6t 
{R — 0.5, r] = 0.02,17 = 1,7 = 1). (b) the simple low-dimensional model 1)5.115.2^ 
(a = 10, /3 = 0.028, 7 = —0.106, S = —0.154). r„ corresponds to the time between two drops and 
Vo is the exit velocity. The two bifurcation diagrams are qualitatively very similar. 



shape after pinch-off. In Fig. ^] by solving H3.6|l numericaUy, we plot the total volume 
V of the drop at the break-up point as well as the ratio of the ejected volume 14 to 
the total volume V as functions of the exit velocity vq. We find that ejected volume 
Ve ~ vqV, so that the volume of the attached drop is a function of the total volume. 
Therefore, we choose the new value of X for the pendant drop keeping dtX the same 
as at the break-up time; this last reinjection process is arbitrary but does not influence 
the results qualitatively. Numerical simulations of the simple model (|5.1I5.2|I plotted in 
Fig. Illb yields a bifurcation diagram for the time between droplets function 
of Vq. Comparing it qualitatively with the same data derived from the hydrodynamic 
model H3.6|l Fig. Illb shows that the transition to chaos occurs in a similar way as the 
flow rate increases, i.e. t„ oscillates as a function of vq followed by a simple period- 
doubling and then by a transition to chaos via a period-doubling bifurcation. When 
the dripping is chaotic we can also notice several qualitative similarities like a period-3 
dripping between two chaotic regime, a reverse period-doubling bifurcation ending in a 
period-2 dripping regime with a large difference in the interval of time between drops, 
etc. The surprising similarities between the bifurcation diagrams suggest that simple 
rationally derived models can capture much of the qualitative dynamics of dripping. 

The model (j5.1l5.2|l is robust with the choice of the functional form a, b, rj and respect 
to the reinjection process emphasizing that the qualitative features of the dripping faucet 
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are due to primarily the Saddle-node bifurcation and the damped oscillations leading to 
it. We note that although our model (|5.1I5.2|1 may be reminiscent of the widely studied 
family of mass-spring models ( see iSKaw 1984|) as well as | |Kiyono et al 1999| ) and refer- 
ences therein) there is a qualitative difference. Unlike in all previous models, we account 
for the all important Saddle-Node bifurcation represented by the nonlinear term in the 
equation (|5.1|l and thus naturally account for the onset of the necking process. 

An even simpler model which ignores the temporal dynamics and focuses exclusively 
on the interval between drop emissions leads to a return map 



dtU ^ {iuj-X)U, 
dZ = e + 



(5.3) 



Here U = X + iY \& the amplitude of the oscillatory damped mode with the associated 
eigenvalue icj — A (A > 0) and e characterises the flow rate and is associated with the 
growing mode Z . As shown in Fig. 112|l the return map is defined using a parallelopiped 
of length (A, A, B) in the phase-space {X, Y, Z) centred at the saddle-node point (Fig. 
p2b )). We first contract the map from the plane Y = before the Saddle-Node area, 
to the plane Z = after the Saddle-Node area. A point {Xi,A,Zi) is mapped into 
(X,+i,r,+i,S) (see Fig. ^)) via 



Xi+i = {XiCoa{ujTi) ~ Asm{ujTi))e 
Yi+i = {XiCOs{LJTi) + Asin{ujTi))e 



where Ti is the return time. 



arctan(i?/-i/e) — arctan(Zi/Y^) 
n = (5.5) 

In order to complete the construction of the dynamical model, we again need a global 
reinjection process that resets the system and replaces the dynamics of pinch-off. The 
simplest way to model the reinjection flow is via a rigid rotation, as for instance 

5'+^^^'+^ (5.6) 
Using (|5.4|) and H5.6|) , the Poincare map which models the process is then given by 

arctan(i?/-ye) — arctan(Zj/\/e) 
= r 

= {XiCOs{ujTi) — Asm{ujTi))e~^'^' ^ ' ' 

Zi+i = (Xi cos(a;ri) + Asin(wTj))e~^'^' 

A numerical simulation of this return map shown in Fig. leads to a bifurcation 

diagram very similar to that of the hydrodynamical model shown in Fig. pib ) and the 
ODE modes (Hit) . 



6. A mechanical analogue of the dripping faucet 

A natural mechanical analogue of the dripping faucet arises from the analysis presented 
in the previous sections and is shown in Fig. 113b . where a frictional particle sits inside a 
corrugated cogwheel rotating about an axis perpendicular to gravity. The particle may 
start at the bottom of one of the local minima, but as the cogwheel rotates, this stable 
static solution eventually disappears via a saddle-node bifurcation. The particle then 
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Figure 12. Bifurcation diagram for the simple return map 15.311 . r„ corresponds to the time 
between two emitted drops and e to the flow rate (Here A — 1.4, B = 1, X — 2.6 and uj = 20). 
Again, we see the similarity to the bifurcation diagram obtained from the complete model (Fig. 
\llh ) as well to that obtained from the ODE model (Fig. Illb'l. 



escapes and performs damped oscillations about a new minimum similar to the previous 
one, and the scenario is repeated. We see immediately that the rotation rate plays the role 
of the flow rate in the dripping faucet, while the periodic potential which is locally cubic 
characterizes the capillary-gravity potential for the drop. The dripping then corresponds 
to the large excursion of the particle from one minimum to another, accompanied by 
inertial oscillations and damping. As the rotation rate becomes similar to the frequency 
of the particle about its local minima, we can expect complex dynamics here just as in 
the dripping faucet. 

To quantify our mechanical analogue, we use an arc-length coordinate S{t) for the 
position of the particle. Then its equation of motion, in dimensionless form, is 



dttS + lydtS = -dsV (6.1) 

Here the vertical position of the particle y{S{t),t) is determined by the instantaneous 
position of the cogwheel which is given as 

x{s,t) = {1 + a cos k {9 (s) + fit)) cos 6 (s) , , 

y(s,i) = {l + acoski0{s) +nt))sm0{s) ^ ' 

As usual, here the relation between the cartesian coordinates and arc-length coordinates 
is given by ds = {dx^ + dy'^y/^ which leads to 

ds = y^l + 2a cos k{0{s) + nt) + cos^ k{9{s) + nt) + a^k^ sin^ k{d{s) + nt)de (6.3) 

In (|6.2|l . n is the angular velocity of the cogwheel, k the wave number of the modulation 
of the cogwheel a the amplitude of the cogs, and is the dimensionless friction parameter. 
When = the particle sits in a stable static position in between two teeth. An unstable 
solution close to the tip of the neighboring tooth separates the particle from another static 
stable solution, and a large enough perturbation will enable the particle to go to this 
solution just as a static drop with a volume less than Vc hanging from a faucet with the 
flow rate q = will drip if perturbed strongly. When $7^0 but is still small, the particle 
moves from one local minimum to the next, accompanied by damped oscillations at a 
typical frequency {dsuY^^, which varies with the location of the particle as for the drop. 
This periodic transition between teeth eventually becomes chaotic when fl becomes large 
enough via a sequence of bifurcations as shown in Fig. 113b. We see that this bifurcation 
diagram is qualitatively similar to those seen in the dripping faucet shown in Fig. and 
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Figure 13. (a) A mechanical model for the dripping faucet consists of a damped particle moving 
inside a cogwheel rotating with an angular velocity fl. The teeth, i.e. the local minima correspond 
to the local minima associated with static drops, and the angular velocity of rotation is similar 
to the flow rate, (b) The bifurcation diagram for the particle in a cogwheel showing r„ the time 
for the particle to move from one tooth to another as function of f2 obtained by solving the 
equation 16.116.21 with k = 13,^ = 0.25, a = 0.01. Again, note the similarity of this bifurcation 
diagram to that shown in Fig. I11I12I 

the hierarchy of models considered earher, i.e. the ODE model bifurcation diagram in 
Fig. 111! and the return map bifurcation diagram shown in Fig. 1121 



7. Discussion 

Based on the study of the stability of a pendant drop and numerical simulations of a 
lubrication- type model for the hydrodynamics of a dripping faucet, we have constructed 
a hierarchy of simple models that qualitatively account for the various experimentally 
observed behaviors of a dripping faucet. Following a numerical simulation of the lubrica- 
tion equations that is used to obtain a bifurcation diagram for the dynamics of dripping, 
we first show that a simple ODE model based on the dynamics of just two modes is suf- 
ficient to capture the qualitative features of the bifurcation diagram. This led us to even 
simpler model of a discrete map which also shows the same qualitative behaviour. Finally, 
we propose a new mechanical analogue for the dripping faucet, one whose behaviour is 
similar to that of the other models. 

All our models differ qualitatively from the many previous models proposed for drip- 
ping over the last two decades, by focusing on the two key elements that govern the 
dynamics of drop formation and ejection : a) the stability analysis of the static drop 
close to its critical volume, b) the reinjection of the system back to the neighborhood 
of a stable static solution. These models are corroborated by a POD of the numerical 
solution determined by solving (|3.t)|l and the spectral analysis of the static solution close 
to the instability threshold. From the perspective of dynamical systems, our models in- 
corporate a coupled system that invokes two classical normal forms associated with the 
Saddle-Node bifurcation and the Shilnikov mechanism. As we have discussed, the former 
corresponds to a loss of a stable "Node" solution via the "collision" with an unstable 
"Saddle" in the context of the linear stability of a drop with V ^ Vc- The Shilnikov 
mechanism for chaos is associated with the presence of damped oscillations and a global 
reinjection IjGuckenheimer and Holmes 19831 lArneodo et al 19 851 which are physically 
manifest in the damped capillary oscillations and the pinch-off process which effectively 
resets the system. 
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Figure 14. Location of the point of drop detachment Zcut as a function of the exit velocity vo- 
We see that as i;o increases, the shape of remaining drop which is originally close to that of the 
stable solution with no waist, becomes closer to that of the first unstable stationary solution 
(the drop with one waist, also shown in Fig. 1). Our simple models are clearly not valid in 
this regime, but by including this mode can be made better. As Vf increases even further, the 
remaining drop becomes long and slender, thus marking the transition from the dripping to the 
jetting regime as seen above. 

We conclude by pointing out the limitations and extensions of our simple hydro- 
dynamical models for dripping. They are valid for small and intermediate flow rates 
when the dynamics of dripping is strongly influenced by the static solution, i..e the 
remainder of the drop after ejection is still similar to the static solution. This notion 
allows us to characterize and explain the transition from dripping to jetting, observed by 
IClanct and Lasheras 1999. Ambravanesw aran ct al 2 004 and others, in a rational way. 
Our simulations of the lubrication equations in §3 showed that the necking time t„ is 
independent of the exit velocity. During this time, the point where the drop eventually 
detaches travels a distance Id ~ voTn (see Clanet and Lasheras 1999). As the flow rate 
increases the length of the remaining drop after ejection increases and so the shape is 
no longer close to one of the stable stationary solutions shown in Fig. ^ At this stage 
the dynamics of the dripping is not influenced by the stable stationary solution, an as- 
sumption which is the main premise of all our models. This defines the transition from 
dripping to jetting; in particular, the influence of the faucet on the dynamics is negligible 
in this regime as the stable static solution ceases to influence the dynamical behavior of 
the system. The long slender fluid filament in the jetting regime also eventually breaks 
into drops through the Savart-Plateau-Rayleigh instability, a regime that is qualitatively 
different from the one we have treated here. 

More generally, our approach should be applicable to a variety of systems where we 
see two relatively common ingredients: the loss of a static stable solution via a saddle- 
node bifurcation and the dynamical return to a nearby state via a damped oscillatory 
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mode. Examples that come to mind readily include the chaotic nucleation of plumes in 
convection, frictional oscillations at a rough interface etc. 
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